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ABSTRACT 


The dynamics of phase locked loops have been studied by many 
investigators for the cases of either a pure sinusoidal signal or 
for a pseudo noise (PN) sequence. This dissertation presents an 
analysis of a phase locked system consisting of a sinusoidal carrier 
modulated with a PN sequence. The effect of the modulation index of 
the phase shift modulation and the design of a loop filter for coherent 
detection of a PN sequence are analysed. Additionally, the error 
probability for each digit, when a likelihood detector is used in 
front of the delay lock discriminator for a PN sequence, is also 


converted to the delay error of the phase lock loop. 
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TABLE OF SYMBOLS 


s(t) = an mesequence with the number of digits Ме?" -1 ina period, 
and the duration of a digit d, amplitude +1 


Ө = the modulation index of the phase shift modulation 
We = the angular frequency of the carrier 
A = amplitude of the incoming signal 


М2 Асо (их++ 05(%)) = incoming signal 


A, = expected value of the signal amplitude at the demodulator 
output 


G,(W)=the power spectral density of s(t) 

n(t)= band limited white Gaussian noise 

n,(t)=an in-phase component of n(t) with respect to the carrier 
n (t)=a quadrature component of n(t) with respect to the carrier 
№(0) = the Fourier transform of n(t) 


- one sided power spectral density of white Gaussian noise 
(watts/ cycle /second) 


F, (w)=the loop filter of the demodulator 


F.(w)=the loop filter of the delay lock discriminator 


В. = band-width parameter of F,(W), (radian/sec) 
Wo = band-width parameter of F,(w), (radian/sec) 
R,(T)=autocorrelation function of s(t) 

R,,(t)=crosscorrelation function of s(t) and n(t) 


E(xJ = the expected value of x 


6;^ = the variance of n(t) 





CHAPTER I 


INTRODUCTION 


A correlation type detector or a matched filter is known to be 
optimum for detecting a signal in additive noise under a wide variety 
of conditions. The noise power is proportional to the effective band- 
width for the case of white Gaussian noise. We can improve the S/N ratio 
by narrowing the effective bandwidth if it does not restrict the infor- 
mation contained in the signal. 


s(t)s(t+ C) 


5(%) + п(+) 3 n(t)s(t * C) y(t) 









Multiplier 





s(t+T) 
Fig. 1-1 A Correlation Type Detector 

Consider the detector in Fig. l-l. s(t) is the incoming signal and 
n(t) is additive white Gaussian noise. If we know the signal except for 
its phase, that is, the phase of the signal is the only information 
contained in the signal, we can generate s(t+T) to multiply with the 
incoming signal. The low pass filter works like an averager in some sense 
as described later in section 4-3, and 

y(t,, Oe eost + Dated cose +T)dt. 
If s(t) and VR assumed to ы. Иса becomes a function 
of only (for long T, and y(t) becomes approximately the sum of the auto- 
correlation function of s(t) "R4(C)" and the crosscorrelation function 
of s(t) and n(t) "RD". 

у(0) = ВОНА, (O 

If we can find a signal such that its autocorrelation function has 
a distinguished difference between the values of T=0 and T#0, we will 
be able to utilize it for phase detection. 
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A sine function {2 A cos Wet has the autocorrelation function of 
К. (С) + А совидт, and the phase locked loop for a sinusoidal function, 
like that shown in Fig. 1-2, has been studied by many investigators. 
Among them, R, Jaffe and E, Rechtin showed an optimum filter design for 
a linearly approximated loop ET A. J. Viterbi showed the transient 
response for 9! having a second order polynomial in the nonlinear 
ѕуѕіет, [10] W. E. Roland defined the threshold criterion by the minimum 
of the cycles of slip anc showed its utility through experimental results 
in his thesis for the degree of H. S. (2 
V2Acos(w.t + Ө)+п(+) ACk, sin(0,-8)* J2Ck,n(t)sin(wct t &) [Low Fass 

Filter 

V2 Csin(wtt %0) 


Voltage Controlled 
Oscillator 


Fig. 1-2 Phase Locked Loop for Sinusoid 









Ihe binary signal has only two possibilities for each digit, ONE or 
ZERO, or -1 ог 31, but sequences of binary digits show many interesting 
and complex characteristics. Some sequences represent codes for which 
maximum likelihood detection is possible. They are classified by the 
crosscorrelation function in each code set; orthogonal, biorthogonal and 
simplex. [2] 

There are also some special group of binary sequences especially 
adapted for phase detection; the PN sequence, especially the maximum- 
length linear shift register sequence (m-sequence), will be explained a 
little more precisely іп Chapter 2. R. H. Barker defined "Group Synchro- 
nization" so as to pin-point effectively the origin in time for a binary 
sequence, and showed a matched filter, and also the characteristics and 
some patterns of this sequence. [1] 


The complementary binary sequences studied by S. Jauregui, Jr. used 
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the same kind of matched filters and when used in a differential pair 
these eliminate small ambiguities in the output. (6] Both sequences can 
be used in the form of discrete pulse sequences. 

PN sequences are adequate for the use of both continuously and 
periodically incoming signals. J. J. Spilker proposed a delay lock 815- 
criminator shown in Fig. 1-3, using the combination of two autocorrelation 
functions having two digit interval difference "7. (8] ‘J. C. Springett 
showed a combination of clock and PN sequence to produce a controlling 
voltage proportional to the phase difference. [9] М. Е. Easterling 
introduced a tracking system using a combination of clock and several codes 


and produced synchonization, code by code, separately. [3] In this system, 





Loop Filter 


Voltage Controlled 
Clock 


Linear Feedback Shift 
Register 


Fig. 1-3 Delay Lock Discriminator for Binary Signal 

only the clock has a phase locked loop. The synchronization of the other 
codes is evaluated by a correlation indicator. For the system which needs 
long distance tracking and does not require quick reaction for the variation 
of the targets dynamics, this will be all right, but if we need automatic 
synchronization and very rapid acquisition, this will not be adequate. 

Since the second one has many weak ambiguous points between distin- 
guished phase lock points, we prefer Spilker's delay lock discriminator 
rather than Springett's. 

Ihe periodic sequence is preferable for phase lock loops to the non- 
periodic because of the possibility of continuous control. 


This paper treats the phase lock system for a carrier phase shift 
11 


modulated by а PN sequence, as shown in Fig. 1-4. s(t) is a maximum- 
length linear shift register sequence with the number of digits in a 

period M-2'-1, and the duration of a digit d. n(t) is white Gaussian 
noise and has the following narrow band approximation 


n(t)- n;(t)cosut * n4(t)sinut 


/2Ac6s(wt+02(t)) AoAlt) 





Fig. 1-4 A Phase Lock System for a Signal Modulated 
by a PN Sequence 


Many papers have analyz2d each separate stage presented so far, but 
not the combined system. There are two problems to be solved for this 
system; what is the effect оѓ Ө оп the phase lock performance of this system, 
and if Е (го) and F,(w) have the same type of low pass filter as that derived 
by Jaffe and Rechtin, how can we determine adequate bandwidths for these 
filters? determines the ratio of the I component and Q component of 
incoming signals, and also the necessary bandwidth for the input of this 
system. If O is large, the output modulating signal will be large, but due 
to the large phase error of the demodulating stage, the effective amplitude 
will decrease. The loop bandwidth of each stage is related to the tran- 
sient response for tracking and the required reduction in the noise output. 

This report will summarize those relations determined by other authors 
for each stage, and will further derive the threshold S/N ratio and the 
performance of this system. Chapter 2 considers the characteristics of 
m-sequences, Chapter 3 shows the demodulation of an incoming modulated 


Е 


carrier and the output of this stage. 

The delay lock discriminator for an m-sequence will be analyzed in 
Chapter 4 with the effect of a limiter in front of this stage included 
in the performance. 

Chapter 5 summarises the system results and illustrates by an example 
the conclusion of the next chapter. This illustration will not show any 
particularly interesting value for@or other parameters. It is simply an 
example of the evaluation of system parameters and of design decision, 
but it is hoped that this will give a rough idea of the overall performance 
of such systems and will give some understanding of the role played by each 
parameter in the design of a system employing PN modulation of a sinusoidal 


carrier. 
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CHAPTER II 


MAXIMUM-LENGTH LINEAR SHIFT REGISTER SEQUENCE 


Ihe »utocorvelation function of a binary sequence of length M, 
3(1)={а,,а,,...,ам|, ућеге а;==1, 15 defined as 

3 енд: auj itj : modulo M 
TO) = 

Ја б)=аф  (Ci«j«n) 
then s(i) is called a twoelevel sequence, and if, moreover, WG- i 
for j#0 (modulo M), the sequence is sometimes called a pseudo-noise 
(PN) sequence, Exactly, a PN sequence has to have the following three 
"randomness properties". [4] 
R=1 (the balance property) In each period of the sequence the number of 
ONE's differs from the number of ZERO's by at most 1, 
R-2 (the run property) Among the runs of ONE's and of ZERO's in each 
period, one-half the runs of each kind are of length one, one fourth of 
each kind are of length two, one eighth are of length three, and so on 
as these fractions give meaningful numbers of runs. 
R-3 (the correlation property) If a period of the sequence is compared, 
term by tern, with any cyclic shift of itself, the number of agreements 
differs from the number of disagreements by at most 1. 


There are four types of PN sequences, corresponding to the length 


м. (2] 

1) M=2” -1 (m-sequences) 

2) M=4t-1 is prime (quadratic residue or Legendre 
sequence) 

3) M=4t-1=4x%27 is prime (Hall sequence) 


4) M=p(p+2) where both p, p+2 are prime (twin prime sequence) 
14 


Among them the maximum-length linear shift register sequence (n-sequence) 
has the delay-and-add property and is generated by a linear shift register. 
This is also convenient to generate in the phase lock loop, therefore in 


this report the m-sequence is used as a FN sequence. 


clock 


ИШҮ: — Out 
fox не, сх 


gc m, уч о пу S, Оці 





Fig. 2-1 n Stage Linear Shift Register 
In the n stage linear shift register shown in Fig. 2-1, the output 
sequence is always periodic with a period p<M=2 -l, and when period is 
M, that sequence is called a maximum-length linear shift register sequence. 
This m-sequence has the delay-and-add property, that is, if the sequence 
is added term-by-term modulo 2 to a cyclic shift of itself, the resulting 
sequence is another cyclic shift of the original sequence. This property 
is possessed only by the binary sequence (а га, „а (а;=1 ог 0) which 
satisfies a linear recursion relation 
n 
ак- до ауы for allk 
k-i : (modulo M) 
с;= 1 ог 0 
summation is modulo 2 
and M=2" 2-1. (2) m-sequence satisfies this. 
Consider an m-sequence s(i)zi(a,, a,,...,a4] - The autocorrelation 
function of this sequence f,(k) is 
5,00 aaa, Ө, sum of modulo 2 
15 


Using the delay-and-add property, 
M 
5, (к) = а, — {CF of ONE's)x(~1)+(# of ZERO's)x(1)} 
ізі 


a оон IL 
= 1 (281 х (-1)+ —ESX(r1)/ 
= = t for k£O (modulo M) 
рші | a 
7, (0) eic -42171 


“ 
£ 


(If we т. ONE and ZERO with -l and +l, we can use multiple instead 
of Q). For the wave form s(t) which takes the value of +1 or -l according 
to s(i) at time (i-l)d<t<id, the autocorrelation function H (T) is, 

? (0) = здае T=Md 
Кок 120 ог ка (k: integer) 

4 2 

$(0)— za ее! dt= 1 

Ӯ, (ка)= m Stt)s(ttka)at = lat > a; ask = - 

If Tekdt+ and O<Ö'<d, then first (d-9) of each digit of s(t) 
corresponds with the last part (d-0) of s(t-c kd), and the remaining Ó of 
s(t) corresponds with the first part Фот s(t+(k+1)d). 


A M 
£4 — arzt а аһ) 
= 1-4 j 


M d M 
%ка «4)- Н [i AES, Gi Gir +), dÈ Z Ai Aeri | 
al 
м 
P, (T) is an even function and we can get the same result for any negative 


value of T, and the autocorrelation function of s(t) has the shape shown 


in Fig. 2-2. 
1 





Fig. 2-2  Autocorrelation Function of m-sequence 
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The power spectral density of s(t), G4(») is derived by using the 
Wiener-Khintchine relation. 
= hi. 
С. (и) << (2 (<) е 


= ны эрш, Hw 211) ++, Sw) 





| || | ава @ме ај 


=» 2! а) 
A = 


Fig. 2-3 Power Spectral Density G,(w) of m-sequence 

When we consider the necessary bandwidth of the band pass filter for 
the modulated input /2Acos(u,t+6s(t)), where O«X0«Z. , we can analyze s(t) 
into frequency conponents. If we pick up one frequency, then that 
frequency component will be 

гені ве. ззп (218. + +Фф,) 
Since M is a "^ nunber, Ма Е «а and the amplitude of each frequency 
component is supposed to be very small. A phase modulated carrier Wwe by 
a modulating frequency Wm, with small modulation index «1, is expressed 
as 

cos(w.ttßsinwt)> cos u.t- Bsinwstsinw.t 
This has frequency components of w,, ёш. Тһе frequency components of 
the signal (2Acos(w.t*tO0s(t)) are derived in the same way as for small 
modulation index «J, therefore in order to preserve almost all power 
of this incoming signal we can pass only the frequency between wti 
because almost all the power of s(t) is included inside of +=" 
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CHAPTER III 
PHASE LOCKED LOOP FOR CARRIER 


3.1 Coherent Demodulator 

Since s(t)= +1, we can separate the input signal \/2Acos(ucttOs(t)) 
into two parts. 

у 2 Асо (ш;%+08 (%) )=Ј2 Асозб0 coswWt+s(t)Y2Asindsinwt 
If we can lock the phase of tne localTosciliator Ро в cw Ner frequens 
we can get the modulating signal V2Asinds(t). While V2AcosÓcos (t can be 
used for that phase lock, Hence, there will be an optimum value of Ó to 
get the maximum signal to noise ratio in the output. 

Phase jitter in the local oscillator, caused by the output of 
(24s(t)sinÓ, plus input noise causes the amplitude variation of signal 
output in the form of \2As(t)sin(O+ 6.) and acts to increase the noise 
input for the next stage. 

The crosscorrelation function of cos(W,t+9s(t)) and cosukt is 

«er = | соз (а5+05(+))соз (++ Т)дЪ 

-z]li [созбсоз (2 wet + W-T)-s(t)sinOsin(2uet + eT) cosw.tcosd 
-Fs(t)sinÜsinaxt ]dt 
= Зсозвсоз оо у sindsinweT 
= Хеозбсов (01 + $ ) 
Ф == tan (ten 8) 

ф is negligible wnen M is large апа tanÜ is not excessive and we can 
make a phase lock loop like Fig. 3-1. 

As far as this phase-lock loop is concerned, the first low pass 


filter can be neglected and letting the low pass filter be denoted as 
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m ва. 27А 
2Acos(u.ttÜ s(t)) X ACk„s(t)sind 


+n(t 
¥2Csinw,t 


L.P.F. F(w) 


Fig. 3-1 Phase Lock Loop for Coherent Demodulation 


Е (0), не can simplify this loop like Fig. 3-2. 


\2Асоз (ес t tri +Os(t) tn(t) ACK, cos0sin(Ai~ Ao)+ACk,,5(t) sinBeos(Ai- Aco) 






-n(t)kY2Csin(wct + X.) 
-V2Csin(w.t+Ao) 


Fig. 3-2 Phase Lock Loop 
Noise input is assumed to be band limited white Gaussian noise, and 
using narrow band approximation 
n(t)- n;(t)eosu.t t n,(t)sinuxt 
where n;(t) and ng(t) are normally distributed time functions varying 
slowly compared withw,, and av.(ni)s av. (n )= O, the variances 
6,,=9,=N="W М bandwidth (cycles/sec) 
7 one sided power spectral density of 
white Gaussian noise (watts sec/cycle) 
When we look at the vector of this input at the vicinity of the phase 
locked point in Fig. 3-3, 2AsinO and ng are found as the causes of 
phase deviation, and since they are independent, their effects will be 


evaluated in the sum of variances. 


(2A 446 
Yl 
On Fon 
У2 А созд У 
E 
Fig. `- ector Liazram for Phase Noise 
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By approximating 

sin(Ai-A)=Ai-Ao 

cos(Aj-A,)= 1 
in the vicinity of the phase locked point, we can linearize the loop as 
shown in Fig. 3-4. The output of the voltage controlled oscillator is 
the derivative of A, and proportional to the input level. Thus we can 


also modify them in Fig. 3-4. 


n 
hit Ost) + ЕВ ID =s ве 40) 
Ло 


Fig. 3-4 Linearized Phase Lock Loop 


(w 
Лу (о) + 05 (w) +e 











+ 


A,(w) 






Fig. 3-5 Linearized Phase Lock Loop 


Let kOe, =k 
AcosO=A, 
then we can get the loop of Fig. 2-5. If we neglect the effect of 
Gs), the solution for optimum performance is shown in the paper of 
Re Jaffe and E. Rechtin. (5] 
Using that result, F,(w) is defined as F,(w)- вен 2 lo апа 


the closed transfer function H(W) is 


H(w) = AKRO с B, + JZ Bw 
/ + АКЕ (оде — В+ /2 В, № + (ји) 


The network for F,(4) is approximately constructed by RC components. 





ма 


К, й 2 
Е К, = Дея -1)» Rz 
Fron the linearized transfer function H(w), we can also derive the 
error Є (уже Мо) Ле (м) 
Лоу 1- НО»? | 
— re - 0? 
xA ) Bi + 42 B,jw — we (3-1) 
For the input of frequency deviation 4W 


Лю 45 





C (u) с (3-2) 
С(6)= == Ее до An etu Et (20 

220 %<0 (3-3) 

0 2 
Е; = | есе 4t = ZEB (3-4) 
Max €(t) is derived by putting dé = : 

Une s (3-5) 
€t). 2145, (3-6) 





The effective bandwidth of H(w) is 
| Г е | 1+2 (96, 
=== и) = u) 
эт) [н ) dw 2m т). | + (з/б, 
-3s- 1.06B, (cycles/sec) (3-7) 
The phase error due to input noise is normal with zero mean and the 


variance is 





6, = zd ol не dw 
7 = 3278. 
sona Z b Ае 80 (3-8) 


3.2 Noise In the Demodulator Output 

As shown in the previous section, the output near the phase locked 
point is 

Acodsin®tAs(t)sindcosd,+n,(t) (3-9) 
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var(n,(t))= zu 

If we subtract the required signal As(t)sinÜ, then the remaining 
is unnecessary and is considered as noise, which is 

AcosfsinB As (6) in ((1-cost) n, (i) 

a) The distribution оғ Ө, 

The phase error 0, in eq.(3-9) is the sum of the contribution from 
the noise input and the self noise of s(t). 

Since Н(ю) is assumed to work in the linear region and ng is 
normal, the phase error produced by the noise input Ө, is approximately 
distributed normally. 

The distribution of phase error by the self noise б, 15 derived 
as follows, 

Consider a randomness property of an m-sequence, "R-2" (the Run 
Property), "Among the runs of ONE's and of ZERO's in each period, the 
number of runs of length i is 1/2 of the runs of each kind, and n= 
log,(M-1) is the maximum length". When we compare these distributions, 
we can see each digit has the probability of 3 to be +1 and $ to be -1. 
If we take a large number of digits m, then the probability it contains 
k of +1 is binomial distribution. 

m 

P(k of +1 out of m digits)=(7)(>) 

Therefore if we take m digits and count the difference of the number of 
-l from the number of +1, the probability of the difference to be i is 

P(difference = i)= (nsi | [ay 
(m and i should be both even, or odd, but since m is large, for a small 
number of i, we can adjust m as m or mtl according to the number of i 
for approximation). For a large number of m, we can approximate it 
with normal distribution. 


Since a low pass filter is, in some sense, an averager, we can 
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conclude that 8 is also approximately normally distributed with zero 
nean and with variance equal to 0 B. d, (See Appendix I). 

0, and 8, have different kinds of input and they are independent, 
therefore, the probability distribution function of 0,= 8.8 is also 
E non with zero mean and variance 6; - ба + 6, : : 
b) The probability distribution об зіп б, апа (1-соѕ 0,) where Q, has a 
i distribution function. 

Considering the noise output by в, if we want to get good per- 
formance in the next stage, Ü should be desired to be a very small value 
in the variance, if not, this system would not work satisfactorily. 

If fhas a small variance, the contribution of 0, бо АсозбзіпВ, і5 
approximately evaluated mainly by 8, neglecting the other terms of series 
expansion of sind, and for (1-со5 ə) by 6/2. 

' AcosÜsinO,- AQcos0 

Ks (t) sánO(1-cos 6,)= ёд. (t)sinO 
The first one has normal distribution and the second one & is Chi- 
square distribution with degree of freedom 1. 

Looking at 0, which is the output of H(w), we can see that it 
varies very slowly compared with each digit interval d, but since the 
ratio of maximum frequency deviation of the input to W.is the same as the 
ratio of the deviation of digit duration to d, and the filter parameter 
В, is wider than the filter parameter «in the next stage, Os varies more 
rapidly than the fluctuation of phase in the next stage. Therefore the 
effect of 6, to the next stage is considered as if Acos@sin§, were added 
to ns. as random Gaussian noise in the approximation of (Асоз0)0, and since 
AsinO0(1-cos0,) is approximated by X,, at the output of the next stage they 
have the same effect as the decrease of amplitude by a factor like an 


average of many digits. 
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The averaged power of As(t)sinÜcosÓ) is 
K sind |, cos рте те, её 6 929, аб, 
– буг ба: -6%64 
E “= Las aat LES (22 6,-6- 40. ] 
ги! + 

— E + Be ("cos 26,2. Z х Ах] 

= sinoli | (3-13) 
Eg 1 
= > Asind 


(6. is distributed in js distribution, and if we take sum of n 
So, 








or (l-cos ĝo)Asinĝ = 


digits, it becomes X. 
£ 
2) = 1, var. (An) = 2/n 

And as n becomes large, the distribution of the average becomes very 
close toie 

àv. (^oinO cos0,) = (1- = )Asin@ (3-14) 
Comparing this with [440% , for small values of бр 

= К z Co 

(је У: ) = l- 82 - 1- 58» 
The result shows good agreement wilh cu. (3-14). 
c) Correlation of Ó, and nə (t) 

From the argument so fsr, we could derive the approximate output of 


tne discrisinstor, and it is shosn as follows. 


s(1+ o в ) (Asinf) (t) as the signal 
+ (Асоѕ0)0, + n.(t) as noise (3-15) 
What is the distribution of noise as a whole? 


o = 8. + Ө» 
Q, is uncorrelated with C and so with n (t). 
Therefore the correlation of 8, and ng(t) is the correlation of 8 апа 
n (t). 


In a linear system of transfer function H(w), 


p Hw) 


h(t) - Еве те cos Be t for t»0 ( h(t£)e0 for t«0) 
24 


бк = а-та 
Correlation of 6), and ng(t ) is 
го А) < geni, n(t+A) A(t)at 
= lin = Š n(t+2) | п(%-т)һ(о)ага 
= ш (Pac) {in(t+A)n(t-CdtaT 
= С най H a(t +A nlt- T)ataT 
= |һ(©в„(А-©ас= | (т NA-Yaz 


= ge^ соз ҒА (3-16) 
for A=0 
rg4(0) — = (3-17) 


По variances of 6, and nalt) are 
22 Bo e 
ба, = = са n + 02É,d 


Or, = (iru (3-18) 


The normalized correlation coefficient Pis 


P= 3200) 


Oe, She (3-19) 
Letting 6,Асозб = х 
no (t) — У 


Е(х) = Е(у) = 0 


бу = А2с05'0 65 


бу = 2% 
-25Х%4, 
Ра NER: | р ё $ 6; £3) 
Changing variables Xty-u 


y= v 
u-vy 2f(u-v)v = 
p(u)— Fan exp [zn 5 (. ex Sx 6y | = дъ 
= есіте 
ехр(- 2(6; TEE ут бу?) 1 я күрт 2566 yrs < dz 
where z= om teeta veg ie > 
| ae + 2 56; o3 
p(u) = 


RV STEER TS d Г z SZ + реву ері 





О, = 6х + Рбхбу + бу“ 

u= n,(t) 
The output is now simplified as 

| (2+ е 4) (asindb(t)+ ng(t) 
n,(t) is normal with zero mean and variance S 

The power spectral density of this noise is not flat in the bandwidth 
of l/d, because 0, has the approximate bandwidth of B,. Let the power 
spectral density of n,(t) be G,(w) 

си ат | б„(е) аш 

If we divide this by дної, , then we can get the effective band- 
width Bg. The low pass filter used in the next stage has narrower band- 
width than Bọ, then we can use Сев) as the approximation of G,(w) because 
the value of G,(W) inside of that bandwidth W» mainly contributes to the 
output. 

Let the Fourier transform of n (t) be М, (№), then the Fourier transform 
of (hcos 8]8, is М (о) (деозб) (и) 

Ма) H(w> „= A 

(Acos ТЕ |. (Асов0) = O 


aw) = {Inet Н) { == 2.917 
ео 0 


Ше 
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CHAPTER IV 


DELAY LOCK DISCRIMINATOR FOR A BINARY SEQUENCE 


4.1 Tracking and Acquisition Performance of a Delay-lock Discriminator 

This has been solved by J. J. Spilker, [8] and here we would like to 
follow the general results. This section explains the tracking and 
acquisition performance, and the next section shows the effect of noise 
on the performance. 

If we take the difference of the outputs from two crosscorrelators 
which differ two digit intervals in Fig. 4-l, and separate it into an 
averaged value A,kD4(€) and another fluctuating part (self-noise) 
A,kn,(t,€), these are derived as follows. 


Loop Filter 


+ 







Input 
A,s(t+T)+n,(t) Bi 


ks(t+T -d) x 


A 
ks (t+T+d) 





Voltage 








"EH сен 
Clock 
n-stage Feedback gain £c 


Shift Register 
Fig. 4-1 Delay Lock Discriminator 
The autocorrelation function of an m-sequence is already shown in 


Fig. 2-1, and from that figure we can draw р,(Є) as this. 


„rel. 
M 





D. (e)= E[s(t+T)s(trarl)-s(t+T)s(t-arb) 
= | R,(€-d)-R,€+d)} (4-1) 
The self-noise in the output of the multiplier network is 
n,(t,e)= s(trdt)s(t+T)=s(ted+T)s(t+T)=D4(€) (4-2) 
For €-kd, k: integer, using the delay-and-add property, 
a) k=0, (modulo М) 
nyt kd) =s(t+dt?)s(t+i) -s(t-d+T)s(t+T) 
=s(t+T+d)s(ttT)-s(t+T-d+d)s(t+T-d) 
=s(t+T+jd)-s(t+T-d+jd) 
—s(trtjd)-s(t*Tr(j-1)a) (4-3) 
b) k=+1, (modulo M) 


n,(t,kd)- s(t«drTrzd)s(t*T)-s(t-arTra)s(t*T)-D4(*a) 





- [1-s(t«T«ia)- MF! 
= .. 
=e(t+Trid)- ——] (4-4) 
c) k-=other integer 

n g(t, kd) = s(ttT+rd)-s(ttT+qd) (4-5) 

Further discussions of this self-noise, i.e. its power spectral 
density and the output from the loop will be explained in the next section. 
A) The output of the crosscorrelator network is 

kAo[Da(e)*n,(t,&)tn (t)/4;] 

This is filtered by a low pass filter grF2(£) and drives the voltage 
controlled clock generator. This relation is now shown in the following 
Fig. 4=2. Obviously the control part of the crosscorrelator output D4(&€) 
is a function of T-T-€, and T is the integral of the controlled frequency. 


T= 5 kA,gc gr F2 (w) BAG )tn,(w,e)tn КОЗІ 
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A,s(trT)« n, (t) кА Па (Еп + по /А 









Crosscorrelator ge Fa (W) 





ks(t+T) 






Fig. 4-2 Simplified Delay Lock Discriminator 


In the locked-on state, Ed, 





Be) - 
л 
€-T-T 
^ A 
jwT — dkA,g. E e FL Є + Лак пода | (4-6) 


Using the same kind of filter as used in Chapter 3, 
Е. (и) = - де + У2 #16) _ ws we + VE Wop 
< 


We + Joj w 4020 
&= KA,g: Se 
^ % 
T . (w^ -Уг иш) Г E По +_No/As 
oe! (4-7) 
The closed loop function H(w) is 
Е Wo + JZ шь) ш 
P ти + V2 wojw + Wwe (4-8) 
ui id тИ nat no / Ao 
4 O " ТУ ] (4-9) 
A 
The value of T/d without noise is 
^ i 
a HOTT 
i T EM — Š 
ae Tt =z (1 негу) = | НИЕМ з шу 1 (4-10) 


Letting те = at 


ч = 4, 





-Wat 
и и Г жез) 
ще , 
0 д с ае — 
max EG = А, t те (4-13) 


The critical value of We for expected maximum value of a is to be 
decided from nonlinear analysis. If we limit this as 
пах. 6/4 <1 


then alu. < 2.19 (4-14) 
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2 Е 2 42 4 
Ет = | eco dt = u, (4-15) 


B) Acquisition transient 





If we take D,(E) only from the previous equation, 


^ 
ED. до (ше + JZ Wejw) M 
à још + 9, (јо) Da C) qr 


M | R Vz ) wy M 
ei + (= — Буе лети (4-16) 
0 е е A 
Let х= 94 y=  y-x-V4 РА ко фи, = ЈЕ 


о, р та .d 
“О, implies mår SET 


Dale) aor is denoted D(x) 


р(х) =) 2-х пе оно. 
x ENDE 
-2-х -2 < x <-1 
(0 otherwise 
aD «D! (x)= П -1 + x < 1 
-l _ü=2 < x < -1 on 1$ х<$2 
0 otherwise 
then (+£ +#=)(-z)= (i+ /2 Rz )DGO (4-17) 
y+ylg=x+x/g+D(x)+ V2D'(x)x 
where х = АХ. х- etc. 
х + (а + 2 (х))х +0(х)- у+у/е (4-18) 


Ihe slope of phase trajectory Y is 


dx .  DOOr([2bDt:/)x- 3⁄4 -4 
| == = ax T > (4-19) 


ns: 


If we plot this performance on the phase plane X and x, D'(x) and 
D(x) are the function of only x and we can find the dividing line on 
Хж-2, Хшеі, x-1l and x2. Let these regions ое called A, B, C, D, E 


as shown in Fig. 4-3. 
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Fig. 4-3 Regions in Phase Trajectory 
Region A and E 


D' (x) = D(x)=0 


X*2X- У+аў " уж em (4-20) 
2 4 + Й 
ER А иа 


$ does not contain x, but is defined only by x and the input. 


Thus, if x=y at the initial point, then 


y= Vx 

Я ee T У: Ж 

u 5 +1 i A. 

ji 

X=ytX, х,: initial position (4-21) 


If x+y at (-0, letting the initial values of x, x, y, y be xe ж, 0, 
Yo» and substitute them into the general solution of this equation 
- T8 
X=c +o, e ^4 y 


and we obtain 


xzytx,g(x, саа (4-22) 
Region B 
D'(x)=-1 
D(x)= -2-x 
__=2-2tl4 - EX - Ha - Й z 
(e ca Р = (4-23) 
so COn t 
Y+ -J2 J 


For a constant velocity input, i.e. y=0, the point having a constant 
Хе k draws a line. 


(kt ~ 2)% -(х+2+у/2)=0 
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If g is large enough to approximate 1/g- 0, then upper part of the line 


combining (x2 -1, x--0.518), (хе -2, Х 0) converges. This line is 





derived as follows, 


Fig. 4-4 Phase Trajectory in Region B 


for constant {=К 


(k-Y2)x - (x+2)=0 


Да a 2-2 It _ 
а a = М2 –- К > 


= = =< — l 
the slope of the line= k = 





therefore k^ -42 k -1=0 
Мз -- / 
к= vz 





„4. l 254 
о yug cp 


Region C 
D(x)= x 
p'(x)2 1 
I e = = қай >» 
Y-- z * (2 4 ет (4-2%) 
Zu xti tu 
7 сауа +/> | 
Region D 
D(x)= 2-x 
D! (x)= -1 


огох: 
y=- 2 Ecc aod d (4425) 
m -24x + 3⁄2. + 
Mp ga J 
Ihe line which divides the plane into converging and diverging parts is 


derived by the same way as that in region B. 





Fig. ^-5 Phase Trajectory in Region D 

The result obtained by CDC 1604 digital computer Гог $ = СО end 
several initial conditions is shown in Fig. 4-6, (Computer programming, 
see Appendix II). 

From this figure x=2.35, x=0 at the initial condition is the limit 
to converge. It implies that © 2.35. 

Time response of position error and velocity error are shown in 
Fig. 4-7, 8. From these figures the locking tine is approximately 
6/w, = 27 [wo 

The phase trajectory, x vs.tTand X vs,7for g=l0 are shown in Fig. 
4-9, 10, 11. From these figures we can find the position error of 0.2, 
and the limit of y is approximately 2.0. Lock time is also 6/w,+2T7 /шь. 
4.2 Noise and its Effect on the tracking 

In this section, I would first like to discuss the effect of self- 
noise; second the effect of the input Gaussian noise, and finally, the 
effect of the limiter before this stage. 
a) Self-noise 

Self-noise is a function of time and phase error. Continuing from 
the last equation (4-2), power spectral density for the time error Є 
being an integer multiple of d is derived as follows; 

n,(t,kd), where k=0 

The autocorrelation function of this self-noise is the sum of three 


triangular waves as shown. By the Wiener-Khinchin relation the power 
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spectral density is the Fourier 


transform of this wave. 








G „© Е=0) 
ee] 
= 2G, (w) (1-созда) (4-26) 
К ==]. 
Gu, (w,e-td)= G,(w) (4-27) 


Since n,(t,d) is the same as signal. 
k=other integer 
By transforming the autocorrelation function 
G,, (» € kd)" 2G, (o) [1-cos (r-q)uq] (4-28) 
We are interested in the effect of this noise in the linearized region 
-4< e <d, 
For the extreme case E=!d, the power ої self=noise in ihe output of 
H(w) is 
б, = | Gn (и, а У|Н соу | "дио 


= (wa ! ға. wdf, - 
EI Lt 20799) {ЧС О Alw - EL) dw 


=T Ja I+ (9/шо) М 
=d Wod (4-29) 
For the least extreme case, 6-0 


б Са о сина (Мама) (in j ual mn Ta | 
vu — 055 red M ШУ 


b) The mean-square delay error caused by the input noise 
= 2. M 50 Gn, (w) 
ба= rn 4 „Ar o | ае 
2. 
M+! Ir (4-30) 


4,3 The effect of input S 








= 3,084” = 


By the previous discussion we approximated the input noise as 
А ° 2 2 š 2 
Gaussian with zero mean and the variance 6,. The signal is a square 


wave of amplitude +A,. 
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Given a digit which has equal probability of +A or -A with additive 
Gaussian noise, and if we decide the signal is rA according to the 
polarity of input, the error probability for each digit is derived as 
follows: 

When we defined A=] +(1 + eg ^A) AsinO, this was derived as an 
averaged value for some interval, hence if we want to get the probability 
of error for each digit, we can not use A, as the amplitude of the signal. 

The probability of error FP. is 

P.= P(s(t)= +1)P(output is negative/s(t)= +1) 

P(s(t)= -1)P(output is positive/s(t)= -1) (4-31) 
The conditional probability of each term is the same value, because 
additive noise has zero mean and the signal amplitude is the same in 
absolute value. 

P.= 2P(s(t)=+1)P(output is negative/s(t)=+1) (4-32) 

The output of the first discriminator is 

Асоз05116, + Аѕ(+) ѕілдсоѕ0,+ no(t) 

Ө. and и; are -Gaussian, with zero mean and 

бә, = ae meng + O Bed 

ET 
they are not mutually independent, but the normalized correlation 
co-efficient fis 

p =EL Se, On, 

From eq. (5-32) 
P(output is negative/s(t) =+1) 
=P [Asin(@ tr 6) + ng < 0] 
—P[n; 4 -Asin(6* &)] 


= 8: 2 pO, ng , Ni 
= тусут), ыы С к + On )j4n; 4 


й b 2 ! 
= Taree J РЕ Ў) т 8,0090 28 ја, (#23) 
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8 
where erf(x)- = j e 2 


Consider a pulse train of m digits, where each digit has probability 
of error Pg, then if all pulses are l, the average of this pulse train is 
derived as follows. 

The number of error k in m digits distributes binormally, and if m 
is large enough and mPg is not small, then this distribution can be 
approximated by a normal distribution with mean=mP., variance 
пР- (1-Р.). In this case each error causes the decrease of 2/m in the 
average, and this decrease is normally distributed. 


Ihe decrease of the averaged voltage v is, thus 


—— es) (4-34) 


If the pulse train is composed with the same number of +l and -1, 


p(v)= 


i.c. n/2 of positive pulses and m/2 of negative pulses, then an error 
which occurred in positive pulses causes negative error 2/m in the 
average voltage and the same for an error in negative pulses. 

X: positive error voltage caused by errors in negative pulses 

y: negative error voltage caused by errors in positive pulses 
Since the error for each digit is independent, x and y are independent 


because they are the error caused by only negative or positive pulses. 
(x)= ==! енто | 
PUT ee actos ie OP | 20Р УМ) 
EM Й Ро НЫ Udo ШЕ 
p(y)= Уат 2 ТЕ (1- Ре)/т exp} ia" 


Since they are independent, the probability distribution of u=(x+ty) is 
со 
p(u)=) Ресор (u-x)dx 
= ср 4 | 
CREER fee = Lom as 
= О Ау) "exp! РЕПЕР ((x РЕУ-к(А “Ре zy )|dx 


т 


І + ` 
= = = \ ехо зе! (Z^ (4-2) Jdéwhere X= Pe =z 


2 alr ee i- io] 
aa e- TIAA 2(2feG- Pa) | (4-35) 
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2Pe (l-Pe) г 
m n S 


This variance is equal to that in the positive pulse train. 

The effects these digit errors cause in the next binary discriminator, 
will be discussed in two phases, first at the locked point, and this can 
be considered as a phase jitter in the phase-locked position, and second 
at the limit where the phase error becomes td. This may be considered as 
the modification of D,(E) curve. These two effects will show the analysis 
of the system. 

When one digit error occurs in the input, how does this appear in 
the outputs of two crosscorrelators? If the system is locked in phase, 
E=Q, the outputs are also a cyclic permutation of that m-sequence and 
their phases differ by 1 digit. Since the autocorrelation function of 
this sequence is nearly zero for é*td, the errors introduced in the 
outputs of the crosscorrelators do not have any relation to their polarity 
except the fact that they occurred at the same time. About the polarity, 
from the characteristics of the randomness of m-sequence, the probability 
that both errors are positive or negative is # each, and the probability 
the sum is zero is 2. This implies the error probability is decreased 
by a factor of 2, and we can use P&/2 instead of Pg, end the change in 
the average voltage output caused by one error becomes double. 

Thus, we can modify the equation of the noise voltage as follows. 

ра) far РЕЗО ED xus ті гриву) | (4-36) 

For the case of ren one error always 522га as negative, and 
another as positive of negative. This tíme we do not need to modify 


Ре, but we can add two effects. 
e = (асай зе 
ЯР ПЕРЕ ОРЕ» expj 
ж 


p(u)= 2 ( £ Ps (1—PE2 ) 
in 
We have been discussing the effect of noise on the average through 
m digits so far, and now we have to make sure of the meaning of the 


average. 
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The inverse Fourier transform of H(W), h(t) is 
- er 

h(t)= 2 ше cos t 
The response of this filter to a pulse of height 1 and duration d is, 
if we take the response after time d, 

wd tt мо фут 

== e ë J2sin( $ t +2) 
Since the chance of error is very small, it will be improbable that the 
digit in error will be discovered. Therefore we can neglect the response 
of this filter in the digit duration and assume this as an impulse response. 

When we look backward for a time interval of T, standing at an arbi- 
trary time t, if we use a real averaging device, the value at t,is, 


assuming ee are k impulses to occur in T, 

u= а t-t: 48-44 , 
If we use ЕС) as an averager, the output average differs depending upon 
the positions of impulses, but when we think of the interval of each 
impulse, this is approximately distributed exponentially, and after taking 
many intervals T containing k impulses, we find k impulses are distributed 
in equal intervals as an average. 
Thus, the output average is the integration of that response. The time 


a; T should be determined for H(w), and looking at the response 


ud o кы (sin ЈИ we can neglect the effect of the tail after t=" 





ү2 © АТА 
r = r 4 (5-37) 
сан з тез N2sin( 2 we t+T)dt 
Rda wë zu ҚМ 
- dicat У, T (au Bt + cove ep) | at 
4 R We en d 
к. аы - 04 ТЕЗ 
ЖЕЛІН = 
n = та 9/74 (4 39) 


Thus we can approximate the effect of H(w) by the effect of averaging. 
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CHAPTER V 


THRESHOLD FOR THIS SYSTEM 


5.1 Summary of the System Performance 

When we evaluate the performance of a phase lock system, maxinum 
frequency deviation, or in other words, maximum target velocity that a 
tracking system can acquire, and transient time for phase lock, probabili- 
ty of cycles slipping and r.m.s. phase error at the steady state will be 
considered. The problem of the cycles slipping is considered very in- 
portant for this system, because the acquisition time in A, E region in 
Fig. 5-3 can be very large. We can not allow cycle slip in the delay 
lock discriminator. Cycle slip is related to r.m.s. phase error at the 
Steady state. From the experimental results of W. F. Roland, [7] we can 


find the linearity of бо, vs. Bo up to 20-30 degree, and the expected 


number of cycle slip per E is 
6.10 * for Cu 20% 
6.10 ^ for бу; 25° 
2.35110 ° for бе : 307 


If we use equation (3-3), the time for SRI. 0.1 is 222.1. Š 
This is for the initial position error of zero, and we can approximately 
use twice the value of it, = , as an е SES time from the initial 
error of iC. 

Let the expected number of cycles slip per 2"/ be p, and the 
settling time be зви then as an average - 

1, = р, 28 ам (5-1) 
is the number of cycles slip in one period. If one cycle of slip occurs, 


þm ë TTC ; 
"в; іп Md is in transient responce and it can cause longer transient 


period of the next stage. 

For the cycles slip in the next stage, Spilker's paper shows that 
for delay error having 6,<0.30d, a probability to lose the locked-on 
state in the absence of transient errors is negligible from experiments. 

While cycles slip determines the limitation for the values of the 
S/N ratio, 0 and the bandwidth parameter B,, Wo, the maximum frequency 
deviation restricts the lower limit of B, and €). From the phase 
trajectory of the transient response of a phase lock loop for a sinusoid 
by A. Viterbi, the meximum frequency deviation at position error = 0 to 
be acquired is approximately 2.15B, for this phase lock system concerned. 
If we restrict ete] € l in the linearized system, 

4g £2.19 from eq. (3-6) 
ar 5572.15 (5-2) 
For the delay lock discriminator, from Fig. 4-6 and 4-10, 
й «2 
This is &< 2. (553) 
Or, using the result of eq. (4-14), 
22.19 (5-4) 
The transient time for phase lock is decided mainly by B, and «e, 


as follows, 


е = = for the phase lock loop 
- 21. 
[ск о 


These transient interval should be limited by the actual situation. 
5.2 An Example of the Calculation Procedure 

We are going to decide Ө, Bə and wə for a radar using this phase 
lock system. The specification of this radar is: 


Range 500 Miles 


^6 


Carrier frequency 200 Megacycles 
Pulse duration, d 5 Microsecond (2500 feet) 
Number of digits in one period M 2-1200 

Т= Ма = 5 Мі1іѕесопа 


Maximum target velocity,v: 3 Mach, 3000 feet/sec 





Input S/N ratio less than 0 dB 
a) Lower limit of Bo, 
Aw| = < We. = 2222 2T. 200x (05 -2tto = Pooo rad aec. 


Bo? =4000 rad/sec 
fomi bhis Bos tp” = = 1.6 msec 
This will be short enough comparing with a period of Md= 5 msec for 
ae first acquisition end we do not need to increase B,. 

For the delay lock discriminator we put 

T(t) =adt in eq. (4-11). 


but tek” Ds (5-5) 
for We=0.6 Чы = 10 sec 

This is too large for this operation, and we need to make wu, large 
in the limit where the cycle slip does not occur. 


b) The number of cycle slips per & 


1,— p, Бе Mà (5-6) 
6 

B,= 4000 d =5.10 M =1200 

= Ар, 


Once a cycle slip occurs, it takes 1.6 msec to recover and this 
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becomes a loss of about one third of a period. Thus we want to limit 
1, 1а бло , or r.m.s. phase error up to 25°, 
с) ба, = 6, + 62, 

ба = mius l.oí B, = amp b (Сања) 

where K is input S/N ratio 


бе = 8,460“ 


2 65 and бо, vs. Ü for various values of K are shown in Fig. 5-1. 


Ас 


~ б%, /. 
ај = ге “уі 


іп“0 

Under the limitation of 60.19 (=(25°) ), we can plot да vs. 6 
for various K, as shown in Fig. 5-2. 

From Fig. 5-2, we can find Ü which maximizes Аг with respect to 
K, and this is in Fig. 5-3. The maximized A2/A*' vs. K is shown in Fig. 5-4. 
e) The mean sequare delay error б, is the sum of the contribution of 
noise and self-noise. Taking the worst value of self-noise, from eq. (4- 
30) and (4-29), 

Si. 3.08 a TEF] + wsd (5-7) 
In Fig. 5-5, we showed буда уз, К with maxinum value of IE, Tor 
various values оҒ 4). In order to avoid cycle slip, we limited б„/а<0.2, 
and the maximum value of “vs K to satisfy this limitation, and t),., vs. 
K for this ware finally shown in Fig, 5-6. 

From these figures, we can see that this system can not show good 
performance for the S/N ratio less than “15 ав. 
5.3 The effect of a limiter before the delay lock discriminator does not 


improve S/N ratio, because for low s/ N ratio, E approaches to + and phase 


lock control voltage D4(&€) effectively decreases to zero. 
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CHAPTER VI 


CONCLUSION 


Ihe result in Chapter 5 shows that if we use the demodulator in the 
limit to avoid cycle slip, the amplitude of the output signal A, is mainly 
decided by sin O0, and r.m.s. phase error Sə, is also determined mainly by 
the noise input near the upper limit оғ 25-30%. Eq. (5-7) shows that the 
variance of delay error 6, is the sum of input noise and self-noise and 
both are proportional to “b, therefore for better S/N ratio the self-noise 
is dominant. We are interested more in the low S/N ratio region than in 
the high one, and we can also approximate 6. by the contribution of the 
input noise. 

When we want to design this phase lock system, we will find several 
opposing factors. 

1) The phase lock ability over the range of the frequency deviation 
limits the lower side of B,, о. 

2) The cycle slip bounds the upper limit of#and B,, JW. with respect to K. 
3) For a satisfactory acquisition time, we have to get large B,and A 
4) The effect of d is related to the input bandwidth and 4e, d is also 
related to the range resolution of a radar. Small d will offer a good 
resolution, but on the other hand it costs S/N ratio and limits the input 
frequency deviation. 

We tried to solve the problem for a constant S/N ratio in Chapter 5. 
Actually input S/N ratio can not be constant in a system application, and 
we can estimate only the least desirable value or a range of S/N ratio; 
When we use this system for ranging, 6, does not have much meaning in 


itself, but is very important as a cause of cycle slip. Looking through 


E 


these four factors, 2 and 3 are éspecially significant in opposing each 
other, and we will finally decide Bo, andwoby applying a weighting 
factor. The actual designing process must be decided according to each 
case, and we want to point out the fundamental relations between the 


above four factors, for they must be further investigated. 
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APPENDIX I 


SELF-NOISE 


The phase error output by signal s(t) itself is 9S(»)H(w) 
considering the relation shown in Fig. 3-5. 
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2i M > rr 2(2 7 a Bo) рү ЫА: 
= 0 M27 <-т-о |+ (ZUBE \ Zu 4. J M 


nto Bd = 
Since M is very large, 44 is very small and we can approximate this 


relation by an integral form. 
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Moreover, let (2) = 26:4 БЕ + 


2 (28+1) (84 )-#< 


and consider contour integrals along the paths a, r, [^ ог б. 


Since Bd is positive, the contributions of the integral path Í to the 


) В,а'ё 


contour integrals | е(2)а2, and з | є(2)е dz, and [% to бе 
с. с 
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/ 


о 
02 


are all zero. The integral along 


і”, 
the path Y does not converge to zero 
for each integral as r-»0, but the 
sum 
[2 


ва ~ Ба = 
gave - [104 a2 -| ве de 


converges as follows. 
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"ы 47 
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For I,, the sum of residues is io 
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If Bod is small enough to neglect second power, 
Ух o e 0 
г-н [150- £089 -a- $5 


="2 Bd = 1.06В,4 


Таз 6:2 06%,4 
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APPENDIX II 


COMPUTER PROGRAM FOR EQUATION (4-18) 


PROGRAM ACQUIRE 
DIMENSION X(30),XD01(30) ,C(15) 
C(10)=1. 

CALL INTEG1(T,X,XDOT,C) 
IF (X(1)+2.0) 10,11,12 
IF (X(1)+1.0) 11,13,14 
IF (X(1)-1.0) 13,13,15 
IF (X(1)-2.0) 16,16,10 
DX-0.0 

рхрот=0.0 

GO TO 2 

DX- -2,0-X (1) 

GO TO 3 

DX=2.0-X(1) 

DXUOT- -1,0 

GO TO 2 

DX-X (1) 

DXDOT-1.0 

ACCE=C(1) 
SPED=C(1)*T+C(2) 
ERROR=C(3)*SPED+ACCE- (C(3)+1.414213*DXDOT)*X(2)-DX 
XDOT(1)=X(2) 

XDOT (2)=ERROR 

GO TO 1 


END 
60 


END 


For g-infinity 
time; from 0. to 9.0, with interval 0.005 
C(1)=0. 
Cli jay 
C(3)=1/g=0. 
Initial value of X(1)--2.0 
X(2)-y 
For g-10 
time; from 0. to 14.0, with interval 0.01 
C(1)=0 
C(2)=y 
c(3)=1/g=0.1 
Initial value of X(1)=-2.0 


X(2)-y 
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APPENDIX III 


SUBROUTINE INTEG 1 


This subroutine provides a simple but moderately flexible means 
for the solution of ordinary differential equations, employing integra- 
tion by the fourth-order Runge-Kutta numerical method. Input, and both 
print and graph output routines are built in. The step size of 
integration can be changed by appropriate data input and/or programming. 

Title: Runge-Kutta Solution of Ordinary Differential Equations 

with Built-in Input and Output Routines 

CO-OP ID: D2-NPS-INTEGL (FORTRAN 60) 

Programmer: J. R. Ward 


Date: June 1964 
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